Atomistic simulations of the incipient ferroelectric KTaC>3 

A.R. Akbarzadeh 1 , L. Bellaiche 1 , Kevin Leung 2 , Jorge Iniguez 3,4,5 and David Vanderbilt 3 

1 Physics Department, University of Arkansas, 

Fayetteville, AR 72701, Arkansas, USA 
2 Sandia National Laboratories, Mail Stop 1421, 
Albuquerque, New Mexico 87185, USA 
^Department of Physics and Astronomy, Rutgers University, 
Piscataway, New Jersey 08854-8019, USA 
A NIST Center for Neutron Research, 
National Institute of Standards and Technology, 
Gaithersburg, Maryland 20899-8562 and 
^Department of Materials Science and Engineering, 
University of Maryland, College Park, Maryland 20742-2115 
(Dated: February 2, 2008) 

Abstract 

A parameterized effective Hamiltonian approach is used to investigate KTaC>3. We find that 
the experimentally observed anomalous dielectric response of this incipient ferroelectric is well 
reproduced by this approach, once quantum effects are accounted for. Quantum fluctuations 
suppress the paraelectric-to-ferroelectric phase transition; it is unnecessary to introduce defects to 
explain the dielectric behavior. The resulting quantum-induced local structure exhibits off-center 
atomic displacements that display longitudinal, needle-like correlations extending a few lattice 
constants 

PACS numbers: 77.84.-s,78.20.Bh,81.30.Dz 
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Numerous experimental and theoretical studies have been carried out on the perovskite 
KTaC>3 over the last forty years (see, e.g, Refs 

[1 [lil [ll] Q an d references 

therein), making this material one of the most-studied "incipient ferrolectrics." The main 
reason for this interest is that the dielectric constant of KTaC>3 increases continuously with 
decreasing temperature down to ~10K, but then saturates to a plateau at a large value (~ 
4000) at lower temperatures while remaining paraelectric and cubic all the way down to zero 
Kelvin . These anomalous low-temperature features are usually thought to be caused 
by the suppression of a paraelectric-to-ferroelectric phase transition by zero-point quantum 
fluctuations (hence the name "incipient ferroelectric" or "quantum paraelectric" used 
to describe KTaOa and other materials, such as SrTiC>3, exhibiting similar unusual dielectric 
and structural properties). Surprisingly, this generally-accepted picture is apparently not 
supported by various first-principles calculations, using density-functional theory (DFT) 
either in its local-density approximation (LDA)[13| or generalized-gradient approximation 
(GGA) 14,ll5f form, since these simulations all predict that KTaO"3 should be paraelectric 

nnn 

at T=0 even when neglecting zero-point motion |4J, y, lf| • This raises the possibility that 
LDA and GGA are not accurate enough to adequately reproduce the qualitative properties of 
incipient ferroelectrics. An alternate explanation for this discrepancy between first-principles 
calculations and experiments is that the simulations assume a perfect material while real 
samples may contain defects such as oxygen vacancies and Fe +3 ions 0, S, Q| that might 
lead to the observed anomalous properties of KTa03. In fact, the interpretations of various 
experiments still remain controversial as to whether they are attributable to extrinsic 

effects (i.e., defects- induced) or intrinsic off-center atomic displacements. Furthermore, while 
previous studies invoke the existence of ferroelectric microregions inside the macroscopically- 
paraelectric KTa03 system to explain some of its properties jl 11|, there has never been 
any direct determination of the size and shape of these proposed polar regions, to the best 
of our knowledge. For instance, the pioneering work of Ref. made several assumptions in 
their analysis of low-temperature Raman spectra - such as isotropy of these microregions - 
to extract a characteristic size ~ 16 A for these polar regions. 

In this Letter, we use large-scale atomistic simulations to shed light on the aforemen- 
tioned long-standing problems. We report calculations on KTa03 using a parameterized 
effective Hamiltonian approach. Our main findings are that (i) LDA and GGA are indeed 
not accurate enough to reproduce the observed anomalous properties of KTa03, even qual- 
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itatively; (ii) these properties can be understood without the need of introducing defects, 
if quantum fluctuations are present to suppress the paraelectric-to-ferroelectric transition; 
(hi) the low-temperature local structure of KTaC>3 is characterized by off-center atomic dis- 
placements that are longitudinally correlated, in a needle-like (and thus anisotropic) way, 
with a correlation length spanning a few 5-atom unit cells. 

We use the effective Hamiltonian (if e ff) approach developed in Ref. [r| to investigate 
finite-temperature properties of KTaC>3. Within this approach, the total energy E tot is a 
function of three types of local degrees of freedom: (1) the u; (B-site centered) local soft- 
mode amplitude in each i 5-atom cell, describing the local polarization in each cell; (2) 
the v i (A-site centered) inhomogeneous strain variables; and (3) the homogeneous strain 
tensor. E tot contains 18 parameters and 5 different contributions: a local-mode self energy, 
a long-range dipole-dipole interaction, a short-range interaction between local modes, an 
elastic energy, and an interaction between the local modes and strains 16]. This effective 
Hamiltonian approach has been successfully used to model, understand, and design ferro- 




electric perovskites (see Refs. [lia, llj], Ilia Ilia, |20j and references therein). E to t is used in two 



different kinds of Monte-Carlo (MC) simulations: classical Monte Carlo (CMC) |2l| . which 
does not take into account zero-point phonon vibrations, and path-integral quantum Monte 
Carlo (PI-QMC) 19, 22, 2^| which includes purely quantum- mechanical zero-point motion. 
Consequently, comparing the results of these two different Monte-Carlo techniques allows 
a precise determination of quantum effects on macroscopic and microscopic properties of 
perovskites. 12x12x12 KTa03 supercells (corresponding to 8,640 atoms) are used in all 
Monte-Carlo simulations. We typically perform 30,000 MC sweeps to thermalize the system 
and 70,000 more to compute averages, except at low temperatures in PI-QMC where more 
statistics is needed. For example, we use 180,000 MC sweeps for thermalization and 240,000 
sweeps at 3K to accurately predict the dielectric response. (Note that we are not aware of 
any previous work reporting the dielectric response computed using PI-QMC) 

In PI-QMC, each 5-atom cell interacts with its images at neighboring imaginary times 
through a spring-like potential (mimicking the zero-point phonon vibrations), while all the 
5-atom cells interact with each other at the same imaginary time through the internal 
potential associated with E tot . The product TP, where T is the simulated temperature 
and P is the number of imaginary time slices (Trotter number), controls the accuracy of 
the PI-QMC calculation. In all our simulations we use TP=600, which we find leads to 
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sufficiently converged results. Outputs of the PI-QMC simulations thus contain local-modes 
Uj(t), where i indexes the 5-atom unit cells of the studied supercell while the imaginary time 
t ranges between 1 and P. Note that CMC simulations can be thought of as corresponding 
to P = 1, so that they do not yield imaginary-time-dependent outputs. 

Figure 1(a) shows the X33 dielectric susceptibility - where the index 3 refers to the [001] 
pseudo-cubic direction - as predicted by the H e g approach, with all its parameters being 
derived from LDA calculations on small supercells of KTaO^ at its experimental lattice con- 
stant. (Technical details of these LDA calculations are similar to those of Ref. It can 
be clearly seen that CMC calculations yield a %33 that is continuously increasing as the 
temperature is decreasing down to nearly zero Kelvin. Turning on quantum effects leads 
to the appearance of a plateau below ~100K with a value of ~100 for the dielectric con- 
stant. These CMC and PI-QMC simulations both predict a cubic paraelectric ground state. 
A plateau for the dielectric response has indeed been experimentally observed in KTa03 
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but reaching a much higher dielectric constant (~ 4,000) and over a much narrower 
temperature range (i.e., below 10K) than in Fig. 1(a). 

In view of these two discrepancies, we have experimented with making minor adjustments 
in the LDA-fitted parameters in the hope of obtaining better agreement with experimental 
data. We have found that this can be done by adjusting just one of the 18 parameters, 
namely, the parameter denoted k 2 in Ref. [l^| . which describes the harmonic part of the 
local-mode self-energy. (In our model, reducing K2 favors ferroelectricity with respect to 
paraelectricity since it leads to a decrease of the zone-center transverse optical frequency by 
weakening short-range repulsions). Figure lb shows that decreasing this single K2 parameter 
by ~18 % from its LDA value of 0.0866 a.u. (atomic units) leads to reasonable agreement 
between our PI-QMC simulations and measurements, not only for the value of the dielectric 
constant plateau, but also at temperatures above 10K. 

Furthermore, this modified K2 also results in a dramatic difference between the two kinds 
of Monte-Carlo calculations. CMC simulations yield a ferroelectric rhombohedral ground- 
state. The corresponding Curie temperature is around 30K, as evidenced by the peak in 
dielectric response displayed in Fig. 1(b). On the other hand, PI-QMC predicts a paraelectric 
ground state. In other words, quantum effects suppress theparaelectric-to-ferroelectric phase 
transition, which is consistent with the accepted picture [JEl- Figures l(a-b) thus (i) reveal 
that extrinsic defects (such as impurities or vacancies), which have been proposed to be 
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responsible for the anomalous properties of KTaC>3 , are not needed to reproduce 

the experimental behavior of this material; and (ii) strongly suggest that, unlike in strongly 
ferroelectric perovskitesj^], IieI. the LDA is not accurate enough for simulating KTaC>3. 

As for GGA, Tinte et al. p report zone-center optical frequencies in cubic KTaOs that 
are all positive and very close to the LDA values. Consequently, according to Fig l(a-b), we 
can conclude that a GGA effective Hamiltonian would not provide a significant improvement 
over our LDA one, and will also fail in reproducing experimental results. This may make 
KTa03 a useful test-case for the development of new functionals within DFT or other ab- 
initio methods. 

We now analyze the microscopic local structure of KTa03 at low temperature. Figure 2 
depicts the magnitude of the local modes Uj inside each i 5-atom cell versus the angle 
that these local modes make with the pseudo-cubic [100] direction, as obtained from a 
T=3K snapshot among the thermally equilibrated Monte-Carlo configurations using E tot 
with the modified Ki- (The magnitude of the local mode is directly proportional to the 
magnitude of the local polarization, e.g., |u| = 0.006 and 0.026 a.u. correspond to a local 
polarization ~ 0.0583 and 0.253 C/m 2 , respectively). Figure 2(a) displays the CMC results, 
while Fig. 2(b) corresponds to PI-QMC 25|. Comparing Figs. 2(a) and 2(b) reveals how 
quantum effects affect the microscopic structure of KTaOs: the local polarizations go from all 
lying close to the [111] direction (corresponding to an angle ~ 54°) and having a relatively 
large magnitude, to being heavily-scattered in direction and having a much smaller but 
non-zero magnitude. The fact that KTaOs is predicted to exhibit non-zero local dipoles, 
even when quantum fluctuations are accounted for, is consistent with the first-order lines 
observed to app ear in Raman spectra which are forbidden in the ideal cubic perovskite 
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Furthermore, an inspection of Fig. 2(b) does not reveal any obvious polar 
microregions. For instance, our quantum-statistical results do not show the local-mode 
distributions breaking up into clusters centered along (111) directions (i.e., angles of ~ 54 
and/or 125°) as would be expected for such polar microregions. 

To gain further insight into the local structure of KTaOs, we decided to compute an 
additional set of coefficients defined as 

iV i=1 \Ui\ |Uj_|_ r | 

Here /i denotes the x, y, or z Cartesian axis chosen along the [100], [010] or [001] cubic 



directions, respectively. The index i runs over all the N B-sites; ttj ijU and Ui+ r ,n are the fi 
components of the local modes in cell i, and in the cell centered at a distance r from cell 
i, respectively. The case in which the local dipoles all have the same (non-zero) magnitude 
and are all aligned along a given (111) direction yields a value of 1 for 9^(r), for any r and 
for any //. This case corresponds to a ferroelectric rhombohedral state having identical local 
and average structures. On the other hand, the other limiting case — for which neighbors at 
a distance r do not exhibit any correlation between the /i-components of their local modes 
— is associated with a zero value for ^(r). 

Figure 3 depicts 9 x (r) (i.e., \i — x) for r lying in the x-y plane. The results correspond 
to one snapshot of a thermally equilibrated Monte Carlo configuration at T=3K, using the 
if e fj with the modified k 2 . Panels (a) and (b) correspond to CMC and PI-QMC simulations 
respectively. One can see that CMC technique leads to a 9 X (r) close to unity for any r, and 
thus generates a macroscopically- and microscopically-ferroelectric rhombohedral structure, 
as consistent with Fig. 2(a). 

On the other hand, PI-QMC simulations give a more complex behavior for 9 x (r) at low 
temperature. One can see that the x components of the local modes are longitudinally cor- 
related in a needle-like fashion: ^(r) adopts large values only when r is along the [100] 
direction. These values decreasing as the magnitude of r increases. (The same result is 
obtained for all symmetry related cases, e.g., for 9 x (r) in the x-z plane, etc.) Figure 2(b) 
further reveals that ^(r) ~ 0.5 for neighbors at a distance of ±2a (where a ~ 4 A is the cubic 
lattice constant) along the x axis. This is in good agreement with the characteristic size of 
16 A extracted from low-temperature Raman spectra of KTa03 jj. On the other hand, our 
simulations go against the hypothesis of isotropic correlation made in ref . jjj] . The longitudi- 
nal needle-like correlations depicted in Fig. 3(b) have also been predicted to occur in classical 
ferroelectrics just above the paraelectric-to-ferroelectric transition temperature j^. In fact, 
they are pretransitional effects that are probably common to most ferroelectric perovskites, 
the peculiarity of quantum paraelectric KTa03 being that the phase transition does not ac- 
tually occur. Finally, note that these needle-like correlations are consistent with the peculiar 



diffuse X-ray scattering observed in Ref. 
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In summary, we have performed large-scale atomistic simulations to investigate 
the (defect-free) incipient ferroelectric KTa03 system using a parameterized effective- 
Hamiltonian approach. The effect of quantum-mechanical zero-point motion is investigated 
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by comparing the results of classical and path-integral Monte Carlo simulations. We find 
that the fitting of all the H e g parameters within LDA yields a theoretical dielectric constant 
that is in poor quantitative agreement with experiment, strongly suggesting that LDA is in- 
adequate for this material. Results in the literature also indicate that GGA will not improve 
the LDA result. On the other hand, a small modification of a single parameter in H e ^ from 
its LDA value is enough to obtain reasonable agreement between theory and experiment for 
the dielectric constant over a wide temperature range. This modified H e $ leads to the pre- 
dictions that (i) KTa0 3 is ferroelectric classically, but becomes paraelectric once zero-point 
phonon vibrations are included, and (2) the quantum-induced local structure of KTa03 is 
characterized by non-zero local dipoles that have longitudinal, needle-like correlations with 
a correlation length spanning a few unit cells. 
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FIG. 1: X33 dielectric susceptibility of KTaC>3 as a function of temperature T. (a) Results for LDA- 
fitted H e fi parameters, (b) Results for modified set of parameters. Solid circles and stars correspond 
to PI-QMC and CMC results respectively. Dashed and dotted lines represent experimental data 
from Refs. j3| and respectively. Solid line shows the fit of the PI-QMC results by a Barrett 
relation A/[(Ti/2) coth(T x /2T) - T ] Q, with A = 27000, T x = 72K and T = 29K. Note that 
our CMC simulations yield a paraelectric-to-ferroelectric transition around 30K, which provides a 
numerical proof for the concept of classical Curie temperature given to To in the Barrett relation. 

FIG. 2: Magnitude of local modes of KTa03 at T=3K versus the angle that these modes make 
with respect to the [100] pseudo-cubic direction. The modified set of H e g parameters is used. 



FIG. 3: Correlation function 6 x {r) of Eq. (0) for KTaOs plotted in the x-y plane for a 12x12x12 
simulation at T=3K. (a) CMC results; (b) PI-QMC results. Each small square represents one 
lattice B site; the origin lies at the center. The modified set of H e g parameters is used. Note that 
the color scales are different in the two panels. 
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